Local prediction of turning points of oscillating time series 
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For oscillating time series, the prediction is often focused on the turning points. In order to pre- 
dict the turning point magnitudes and times it is proposed to form the state space reconstruction 
only from the turning points and modify the local (nearest neighbor) model accordingly. The model 
on turning points gives optimal prediction at a lower dimensional state space than the optimal lo- 
cal model applied directly on the oscillating time series and is thus computationally more efficient. 
Monte Carlo simulations on different oscillating nonlinear systems showed that it gives better pre- 
dictions of turning points and this is confirmed also for the time series of annual sunspots and total 
OO ' stress in a plastic deformation experiment. 

o 

o . 

fSJ , The prediction of oscillating time series that do not exhibit apparent periodicity has been a long lasting challenge and 
^jf^' the focus of three time series prediction competitions P, 0, Q . Evidence from the winning models of the competitions 
^ , and other prediction studies raises two main points: (a) multi-step ahead prediction requires a long time window and 
consequently a high embedding dimension M, and (b) local prediction models, called also nearest neighbor models, are 
computationally efficient and compete other more complicated black-box models, such as neural networks. Notably, 
they were among the winning entries of the two first competitions. In this work these two points are incorporated in 
the prediction of the turning points of the oscillating time series. The prediction of successive samples is associated 
with intra-oscillation correlations whereas the prediction ofjturning points regards inter-oscillation correlations, which 

Uare more relevant to the underlying oscillating dynamics \^\. Turning point prediction is of great practical interest 
. ' in many applications, such as finance 0, Q. It will be shown below that the prediction of turning points with local 
P5 models can be improved using state space reconstruction solely on the turning points at a lower embedding dimension 
: m. 

Ph . For an oscillating time series of length N , {x{t)}^^-^^, where the observation time is tTs and Ts is the sampling time, 
a sample x{t) is a turning point if it is the minimum or maximum of all samples in the time window [t — p,t -\- p] 
I where the parameter p determines the tolerance for temporal closeness of successive turning points. A small p may 
^ ■ assign turning points for glitches in the case of noisy oscillations, whereas a large p may not detect peaks and troughs 
. of short lasted oscillations. For noisy time series, a small p can still be used in conjunction with filtering, and then 
' the turning points are located on the smoothed time series. Note that the turning point magnitudes are then taken 
from the original time series. Denoting the turning point t/i = x{ti) at time point ti, we derive the time series 
and {tiY^^i of magnitudes and times of the alternating turning points, respectively. Thus two successive samples of 
{2/i}"=i regard an oscillation of {x{t)}^^-^. One may also consider three turning points to include both the start and 
. end of the oscillation. 

' In the reconstruction of the Af -dimensional pseudo-state space from {x{t)}^^-y, the reconstructed points have the 
. . , general form x(i) = [x{t), x{t — ti), . . . , x{t — tm-i)]'- The standard delay embedding suggests the use of a properly 
^ ■ selected fixed delay r, so that the time window length is = tm-i = {M — 1)t. For multi-step prediction, should 
be large enough to account for the mean orbital period of the underlying to time series trajectory and it should cover 
the period of an oscillation or a pattern of oscillations [tI] . The choice of large delays or a fixed large r instead of a large 
' M would not be appropriate as in this case large pieces of information from the oscillation, most importantly the peak 
and trough, may not be represented in x(i). The main idea in the proposed approach is to let Tj, j = 1, . . . , M ~ 1, 
and r„j vary with the target time t, so that the peaks and troughs are selected as components of x(i), resulting in a 
smaller embedding dimension m than M. Moreover, the reconstructed trajectory {x(t)}^^^ is subsampled at times 
This is actually the state space reconstruction of {yi}2=i in an m-dimensional state space. The reconstructed 

point is 

Yi = [yi,yi-i,---,yi-m+i]' 

= [x{ti),x{ti^l), . . . ,x{ti-rn+l)]', 

for i — m, . . . ,n, where the implied lags with regard to x(t) are Tj — ti ~ ti^j, j = 1, . . . , m — 1. 

The compression of {x{t)}^^i to {yi}f^i simplifies the embedding because m < M and no lag parameter is involved, 
at the cost of stripping off the information in the samples between the turning points. This compression is analogue to 
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the reduction of a flow to its Poincare map. Note that formally Poincare maps require that a state space reconstruction 
is made first, whereas in this approach the state space reconstruction is made on the turning points. For low- 
dimensional chaotic systems with sheet-like structure, it has been shown that the local maxima, i.e. every second 
turning point, reproduce the dynamics of the respective Poincare map that has fractal dimension one less than that of 
the flow, e.g. the time series of the third variable of the Lorenz system iSj] (see also j9j for the so-called peak-to-peak 
dynamics). Similarly, we expect that the time series of magnitudes of turning points preserve the original dynamics 
of the flow and form an attractor with a fractal dimension smaller than that for the flow and somehow larger than 
that for the corresponding Poincare map, as estimated through the local maxima. This is demonstrated in Fig. [T] for 
the Mackey Glass delay differential equation with delay A = 30 that regards a correlation dimension v ~ 3.0 jlO| . 
The delay differential equation is solved with a discretization step of O.ls. The local slope of the correlation integral 
estimated on a densely sampled oscillating time series, i.e. Tg = 5 s, = 200000 and roughly 3200 oscillations, does 
not maintain sufficient scaling and the same holds for the respective time series of turning points (see Fig. [IJi and 
c) . The insufficient scaling persists for the oscillating time series even for a larger sampling time (ts = 30 s) giving 6 
times more oscillations for the same N (see Fig. [T|d), whereas scaling at the level v ~ 2.2 is formed from the turning 
point time series of the same iV, as shown in Fig. [TJi. 

(a) (b) 




FIG. 1: Local slope vs base 10 logarithm of distance r for embedding dimensions 1, . . . , 10, as indicated in the panels, and for 
oscillating time series in (a) and (b) and turning point time series in (c) and (d) from the Mackey Glass system with A = 30. (a) 
r, = 5s, iV = 200000, r = 35; (b) = 30s, = 200000, r = 6; (c) r, = 5s, iV = 200000, n = 6380; (d) r, = 5s, iV = 6284867, 
n = 200000. The lag r is selected from the minimum of mutual information. 

The amount of loss of information in the compression of {x(t)}^]^ to {yij^Li and depends on the curvature 

of the upward and downward pattern of the oscillations. Actually, in the case of linear upward and downward 
trends there is no loss of information, as each sample x(ti — k) between two turning points x{ti^i) and x{ti), where 
k G {0, 1, . . . , — ^i-i}, can be expressed in terms of the magnitude and time of the two turning points as 

X{U - k) = X{U) - fci^M^^ . _ k^l^. 
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The delay embedding on {yi}2=i does not account for the times of the turning points and information from 

the samples is lost. However, for prediction purposes, it is operationally tractable to form the reconstructed state 
space from the magnitudes {uiYi^i in order to find neighboring points for the local prediction scheme and then call 
in the times {ti}f^i to estimate the time position that corresponds to the predicted magnitude of the turning point, 
as shown below. Our attempts on simulated chaotic systems with dynamic local regression models making use of 
magnitudes and times to reconstruct the state space showed no improvement in the prediction of turning points. 

The prediction model of choice in this work is the local average mapping (LAM), but other local models can be 
developed in a similar way. For a fixed number of neighbors K and given the turning points up to time t^, the one-step 
ahead prediction of the turning point magnitude yi+i is estimated by the average of the one-step ahead mappings of 
the K nearest neighboring points yi{k), fc = 1, . . . , /C, to the target point 

1 ^ 

TJi+l ^ -^^y^(k) + l■ (2) 

The prediction of the time of the turning point yi+i, is estimated from the average of the corresponding time 
increments of the K neighboring points 

1 ^ 

ii+i = ti + —'^^{tn^k)+i - ti(k)). (3) 
fe=i 

For the iterative prediction of the turning point magnitude at a lead time T, the target point at time t^+i is updated 
as Yi+i = [yi+i, Vi, ■ ■ ■ , yi-m+2]' , and the one-step prediction is done as in ([2|) but for the new set of neighboring points 
of Yi+i- This step is repeated until the prediction of yi+r is reached. The direct prediction scheme is simpler and 
faster as it predicts yi+T directly by the average of the T-step ahead mappings of the neighboring points of yi . The 
iterative and direct multi-step ahead prediction of the times of the turning points is done similarly using ([3|) . We refer 
to this prediction model as extreme magnitude local average map (EMLAM) to stress that the neighboring points in 
the local average map are formed only on the basis of the magnitudes of the local extremes (turning points) . 

We investigate whether we can predict the forthcoming local extremes of an oscillating time series better than 
standard multi-step prediction with LAM. In Fig. [2] an example is shown for multi-step prediction of turning points 
of the fourth variable of the Rossler hyperchaos system The predicted turning points with LAM, are identified 
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FIG. 2: Out-of-sample direct predictions of turning points for the fourth variable of the Rossler hyperchaos system {js = 0.1s) 
with LAM (M = 9, r = 10, iiT = 5) and EMLAM (m = 3, iiT = 5) as given in the legend. The vertical line is at t = 15000 of 
the current turning point set to for clarity. For LAM t is advanced by p = 3 to account for the time the target turning point 
is detected. 

from the multi-step sample predictions in the same way as the turning points are determined on the oscillating time 
series. This time series has a rather linear upward and downward trend, so that the loss of information using only 
the turning point time series is expected to be small. 

The superiority of EMLAM over LAM for this system is confirmed by simulations on 1000 realizations using the 
prediction measure of normalized root mean square error (NRMSE) on the last quarter of each time series. As shown 
in Fig. [31 the prediction with EMLAM is better both for the magnitude and time of the next turning point and this 
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FIG. 3: The average NRMSE of the prediction of next turning point with LAM (direct scheme) and EMLAM {K = 5) from 
1000 time series of length 8192 from the fourth variable of the Rossler hyperchaos system. The error bars denote the standard 
deviation of NRMSE. The in the abscissa is defined as r„ = (M - 1)10 for LAM and = [m - 1)33 for EMLAM, as the 
mean oscillation period is estimated from the power spectrum peak to be 66. In (a) and (b) the data are noise-free and the 
prediction is for the magnitude and time of the turning point, respectively. In (c) and (d) the same predictions are for data 
corrupted with 10% observational white normal noise. In order to detect the turning points a zero-phase filtering of order 13 
is applied to the noisy data. 



holds for noise-free and noisy data. The difference is smaller for the turning point magnitudes of the noisy data. The 
difference of LAM and EMLAM prediction is larger when the iterative scheme is used for predictions with LAM. The 
superiority of EMLAM over LAM prediction persists for different data sizes iV, number of nearest neighbors K , and 
prediction steps T, as shown in Table [J Table U shows the least NRMSE and the corresponding embedding dimension 
for the range of the other factors. These are iiT = 1, 5, 10, T = 1, 2, 3, = 4096, 16384, and are shown in the rows, 
whereas magnitude and time are in the columns. 

EMLAM provides computationally efficient predictions at a small embedding dimension up to to = 4, whereas LAM 
fails, at cases dramatically, to reach the level of prediction of EMLAM for any of the tested embedding dimensions. 
The direct prediction scheme shows less dramatic differences in the performance of the two prediction models. With 
the addition of observational noise the results are qualitatively the same and the differences get smaller for the 
direct prediction scheme (see also Fig. 3) and larger for the iterative prediction scheme. These results are based on 
simulations with 5% and 10% observational noise, not shown here. For the detection of the noisy turning points, 
zero-phase filtering was used with an order adjusted to the amount of noise in the data, in order to smooth out close 
local maxima and minima that apparently do not correspond to real oscillations. 

The same simulations have been applied to other oscillating time series of varying complexity that are not charac- 
terized by linear upward and downward trend, namely the first variable of the Rossler hyperchaos system, the first 
and third variable of the Rossler system [12], and the Mackey-Glass delay differential equation for delay parameter 
17, 30, and 100 [l^. The overall results show that EMLAM gives as good, and at cases better, predictions of turning 
points as the ones obtained by LAM. 
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TABLE I: Summary results of the average NRMSE as for the noise-free case in Fig. |3]but for varying A'^, T and K. For each 
combination of A*', T and K, the M of best prediction with LAM and m of best prediction with EMLAM together with the 
respective NRMSE are given, where M = 3, . . . , 10 (t = 10) and m = 2, . . . , 6. 

The simulations revealed some important features of turning point prediction in favor of EMLAM. In all cases, the 
best predictions with EMLAM were obtained with a small m at the level of the fractal dimension of the underlying 
system, e.g. for the Mackey-Glass system with delay 100 that has a fractal dimension about 7, best results were 
obtained for m at the range from 7 to 10. For LAM, best results could be reached only for large r„ implying very 
large M . Another interesting feature is that for LAM the direct scheme predicts the turning points better than the 
iterative scheme, whereas for EMLAM both schemes give similar predictions. It is also noted that in the noise-free 
case, a smaller sampling time r, improves the accuracy in the detection of the turning points and consequently 
enhances the prediction with EMLAM, whereas it perplexes the selection of the embedding parameters M and r for 
LAM. 

We apply the same prediction setup to the celebrated time series of annual sunspot numbers from year 1700 to 
2006. Sunspot numbers exhibit a rather regular oscillation of about 11 years long with stable trough but varying 
peak that has given rise to debatable prediction results suggesting stochastic, noisy periodic, and chaotic behavior, 
obtained with time series models and other models, such as models of input-output systems [l|, lilQEIill- As 
shown in Fig. [4^, the out-of-sample predictions of LAM give rather periodic oscillations, failing to approximate the 
true peaks, whereas EMLAM matches better the true peaks both in magnitude and time. The difference in LAM 
and EMLAM turning point prediction is rather consistent over different embedding schemes, as shown in Fig. 
regarding the most favorable scenario for LAM. Other state space reconstructions with r > 1 as well as the direct 
scheme gave worse LAM predictions. The summary results in Table |TT] for different K and T, show the superiority of 
EMLAM over LAM, where again best EMLAM predictions are obtained for small m. 

Finally, we compare LAM and EMLAM on the time series of total stress from an experiment of plastic deformation 
that exhibits the Portevin-Le Chdtelier (PLC) effect. Poly-crystal Cu-15%A1 is tensile strained at e = 6.67 • 10^^s~^ 
and T = 125°C, and the total stress is sampled at Tg = 0.2 s; further details can be found in [l,[l3]- The increasing 
trend of total stress was removed and the predictions were done on the last quarters of overlapping segments of 
duration 1200 s and sliding step 300 s. For a range of model specific parameter values EMLAM gave consistently 
better predictions than LAM. In Fig. [5l the turning point magnitude and time predictions with the direct scheme are 
shown for LAM and EMLAM for arbitrary chosen model parameters. EMLAM improves drastically the prediction 
of magnitudes and times of the turning points for most of the epochs and the difference is larger for the iterative 
prediction scheme (not shown here). LAM predictions are essentially at the level of mean value prediction. To the 
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FIG. 4: Iterative prediction of turning points of sunspot numbers with LAM and EMLAM. (a) Out-of-sample iterative prediction 
of the last 70 annual sunspot numbers with LAM (M = 11, t = 1, K = 5) and EMLAM (m = 3, K = 5) as given in the legend. 
The turning points were detected using p — 1. The vertical line denotes the target time, (b) NRMSE of iterative prediction of 
turning point magnitudes at T = 1, 2, 3 with LAM and EMLAM {K — 5), as given in the legend, computed on the last quarter 
of the sunspot time series, where r^, = M — 1 for LAM and r^, = (m — 1)5 for EMLAM to account roughly for the 11 year 
cycle. 
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TABLE II: Summary results of NRMSE for direct prediction of sunspot turning points structured as in Table U where M = 
6, 11, . . . , 46 (r = 1) and m = 2, . . . , 10. 

contrary, EMLAM attains much smaller NRMSE that varies in a way that allows the identification, at some degree, of 
the transition of PLC band types: from type A to type C and then back to type A, as shown with the vertical lines in 
Fig. [S] The actual transitions of PLC bands could only be identified by special equipment following the local strains 
along the specimen (for details see [13 )■ So, beyond improving the LAM prediction, EMLAM prediction can possibly 
be used as a discriminating measure for the PLC band types. This point certainly bears further investigation. 

In conclusion, this work suggests that the analysis of oscillating time series, in particular those exhibiting rather 
linear upward and downward trends, can be improved, and simplified, by restricting the analysis to the turning points. 
It was shown that the information in the time and magnitude of the turning points can be adequate to explain the 
system dynamics. Simulations on a number of chaotic flows and the two real-world examples showed that a local 
average model based only on the turning points can predict turning points equally to, or better than, the standard 
local average model, a result of paramount importance for many applications. There are a number of issues to be 
addressed, such as implementation of other model types and inclusion of turning point time in the model, but it seems 
that the focus on turning points can give a new perspective in the analysis and prediction of oscillating time series. 
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FIG. 5: Direct prediction of the magnitudes of turning points of total stress in consecutive overlapping segments at times given 
in the abscissa. The detection of turning points is done with window of p = 1 and zero-phase filtering of order 13, and the 
model parameters are K — 10, r = 5, AI — 10 and m — 5. The predictions are for one and two time step ahead for LAM and 
EMLAM as given in the legend. The vertical grey lines denote the separation of PLC band types denoted at the bottom, (b) 
The same as in (a) but for the times of the turning points. 
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